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We present an algorithm which, given such a group F, returns a fundamental 
domain and a finite presentation for F with a computable isomorphism. 
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Introduction 



Starting with a suitable quaternion algebra over a number field F with exactly 
one complex place, one can construct discrete subgroups of PSL2(C). These groups, 
called arithmetic Kleinian groups, act properly discontinuously with finite covolume 
on the hyperbolic 3-space. Our main result is: 

Theorem 1. There exists an explicit algorithm fAlaorithm \lS\ with Alaorithm \10\ 

for enumeration) which, given an arithmetic Kleinian group T described by an order 
in a quaternion algebra, returns a fundamental domain and a finite presentation 
for r with a computable isomorphism. 

It has applications both in hyperbolic geometry and number theory. On the 
geometrical side, it provides explicit description of a large family of hyperbolic 3- 
orbifolds, which are still not well understood. The algorithm described here can be 
used to experimentally investigate conjectures about them. On the number theo- 
retical side, the algorithm presented here prepares the ground for computing the 
cohomology of these groups with the action of Hecke operators, which gives a con- 
crete realization of certain automorphic forms |Fra98] . By the Jacquet-Langlands 
correspondence |JL70j . such forms are essentially the same as automorphic forms 
for SL2 /F. They should have attached Galois representations, but the construction 
of these representations in general is still an open problem. Our algorithm could 
also allow for empirical study of these objects. 

This problem has already received some attention. In the analogous Fuchsian 
group case (a subgroup of PSL2(M)), an algorithm may have been known to Klein. 
J. Voight |Voi09| has described and implemented an efficient algorithm exploiting 
reduction theory. In the special case of Bianchi groups, that is when the base field 
is imaginary quadratic and the algebra is split, R.G. Swan |Swa71| has described an 
algorithm, that was implemented by Riley [Ril83| and A. Rahm I RahlO| : D. Yasaki 
[YaslOj has described and implemented another algorithm based on Voronoi' theory. 
C. Corrales, E. Jespers, G. Leal and A. del Rfo |CJLdRd4] have described an 
algorithm for the general Kleinian group case. They implemented it for one division 
algebra with imaginary quadratic base field. Our algorithm is more efficient and our 
implementation more general. We have recently found an unpublished algorithm of 
K. N. Jones and A. W. Reid, mentioned and briefly described in ICFJROll section 
3.1], which also solves the same problem. 

The article is organized as follows. In the first section we recall some basic 
definitions and properties of hyperbolic geometry, quaternion algebras and Kleinian 
groups. In the second section we describe our algorithms: basic procedures to 
work in the hyperbolic 3-space, algorithms for computing a Dirichlet domain and a 
presentation with a computable isomorphism for a cocompact Kleinian group, and 
how to apply these algorithms to arithmetic Kleinian groups. In the third section 
we show some examples produced by our implementation of these algorithms and 
comment on their running time. 



I would like to thank John Voight for proposing me this project and supervising 
my Master's thesis, and Karim Belabas and Andreas Enge for their helpful com- 
ments on earlier versions of this article. Experiments presented in this paper were 
carried out using the PLAFRIM experimental testbed, being developed under the 
Inria PlaFRIM development action with support from LABRI and 1MB and other 
entities: Conseil Regional d'Aquitaine, FeDER, Universite de Bordeaux and CNRS 
(see https : //plaf rim .bordeaux. inria. fr/p . 
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1. Arithmetic Kleinian groups 



Here we recall some basic definitions and properties of hyperbolic geometry, 
quaternion algebras and Kleinian groups. The general reference for this section is 
[MH03] . 

1.1. Hyperbolic geometry. The reader can find more about hyperbolic geometry 
in |Rat06j . The upper half-space is the Riemannian manifold T-L^ = C x R>o with 
Riemannian metric given by 

, P dx^ + dw^ + dt^ 
- ^ 

where {z,t) g 'H^ , z = x + iy and t > 0. For w,w' G H^, d{w,w') is the dis- 
tance between w and w'. The set P"'^(C) is called the sphere at infinity. The upper 
half-space is a model of the hyperbolic 3-space, i.e. the unique connected, sim- 
ply connected Riemannian manifold with constant sectional curvature —1. In this 
space, the volume of the ball of radius r is 7r(sinh(2r) — 2r). 

The group PSL2(C) acts on in the following way. Consider the ring of 
Hamiltonians M = C + Cj with multiplication given by = —1 and jz = zj 
for z € C, and identify with the subset C-|-M>oj C H. Then for an element g = 
(2 ^) e SL2(C) and w G 'H^ the formula 

g ■ w ^ {aw + b){cw + d)^^ = {wc + d)^^ {wa + b) 

defines an action of PSL2(C) on by orientation-preserving isometrics. This 
action is transitive and the stabilizer of the point j £ in PSL2(C) is the sub- 
group PSU2(C). 

The trace of an element of PSL2(C) is defined up to sign, and we have the 
following classification of conjugacy classes in PSL2(C): 

• If tr((7) G C \ [—2, 2], then g has two distinct fixed points in P^C, no fixed 
point in V.^ and stabilizes the geodesic between its fixed points, called its 
axis. The element g is conjugate to ± ^"i) with |A| > 1; it is called 
loxodromic. 

• If iT{g) G (—2,2), then g has two distinct fixed points in P^C, and fixes 
every point in the geodesic between these two fixed points. The element g 
is conjugate to ± (^^'^ ^^^g ^ with G K \ (tt -I- 27rZ); it is called elliptic. 

• If tr((7) = ±2, then g has one fixed point in P^C and no fixed point in Ti,^ . 
It is conjugate to ± ( J } ); it is called parabolic. 

1.2. The unit ball model. In actual computations we are going to work with 
another model of the hyperbolic 3-space. The unit ball B is the open ball of center 
and radius 1 in = C -I- Rj C H, equipped with the Riemannian metric 

2 _ 4(da;2 -|-dy2 + di2) 
~ (1- |w|2)2 

where w — {z,t) G B, z = x + iy and jwp = x^ + y^ + t^ < 1. The sphere at 
infinity dB is the Euclidean sphere of center and radius 1. The distance between 
two points w,w' €z B is given by the explicit formula 

diw, = cosh- (l + 2 (^_|i;;2;(rl'Kp) ) • 
The upper half-space and the unit ball are isometric, the isometry being given by 

(n^ ^B 

[ wi — > {w - - jw) ={l-wj) (w-j), 
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and the corresponding action of an element g = ( ° ^ ) G PSL2(C) on a point w G B 
is given by 

(1) g ■ w ^ {Aw + B){Cw + D)-^ 

where 

A = a + d+{b- c)j, B = b + c+{a- d)j, 
C ^ c + b+ {d- d)j, D = d + a+ {c- b)j. 

In the unit ball model, the geodesic planes are the intersections with B of Eu- 
clidean spheres and Euclidean planes orthogonal to the sphere at infinity, and the 
geodesies are the intersections with B of Euclidean circles and Euclidean straight 
lines orthogonal to the sphere at infinity. A half-space is an open connected subset 
of B with boundary consisting of a geodesic plane. A convex polyhedron is the in- 
tersection of a set of half-spaces, such that the corresponding set of geodesic planes 
is locally finite. 

1.3. The Lobachevsky function and volumes of tetrahedra. We are going 
to compute hyperbolic volumes, and for this the main tool is going to be the 
Lobachevsky function, which we define here. The integral 

In 12 sinwl du 







converges for 6 E M. \ vrZ and admits a continuous extension to R, which is odd and 
periodic with period tt. This extension is called the Lobachevsky function C{9). The 
Lobachevsky function admits a power series expansion, converging for e [— tt, tt]: 

With this function one can derive a formula for the volume of a certain standard 
tetrahedron. We will use it to compute the volume of convex polyhedra. 

Proposition 2. Let T be the tetrahedron in Ti.^ with one vertex at go and the 
other vertices A, B, C on the unit hemisphere such that they project vertically onto 
A',B',C' in C with A' = to form a Euclidean triangle, with angles ^ at B' and 
a at A' , and such that the angle along BC is 7. Then the volume ofT is finite and 
given by 
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'TT 

,2 ^ 

Proof. This formula can be found in |MR031 paragraph 1.7]. □ 



Vol(r) = - C{a + -t) + C{a - 7) + 2/: 



1.4. Kleinian groups, Dirichlet domains and exterior domains. A sub- 
group F of PSL2(C) is a Kleinian group if it acts discontinuously on , or equiv- 
alently if it is a discrete subgroup of PSL2(C). A fundamental domain for V is an 
open subset F of 'H'^ such that 

(i) U^er7:^ = H^ 

(ii) For aU 7 e F \ {1}, J" n 7J" = 0; 
(ui) Vol((9J") = 

where Vol is the Riemannian volume on Ti.^ . To compute a fundamental domain 
for a Kleinian group F, we are going to use the standard construction of Dirichlet 
domains. Let p E B he & point with trivial stabilizer in F. Then the Dirichlet 
domain centered at p 



-Dp(F) = {x G S I for all 7 G F \ {1}, d(x,p) < d(7a;,p)} 
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is a convex fundamental polyhedron for F. If F has finite covolume, then Dp(T) fas 
finitely many faces. A Kleinian group F is geometrically finite if one (equivalently, 
every) Dirichlet domain for F has finitely many faces. 

Note that since F acts properly discontinuously on B, almost every point in B 
has a trivial stabilizer in F. In the unit ball model, the Dirichlet domain centered 
at has a simple description. Suppose g e SL2(C) does not fix in B. Let 

• 1(g) ^{weB\ d{w, 0) d{gw, 0)}; 

• Ext(g) ^{w<eB\ d(w,0) < d(gw,0)}; 

• Int(g) = {w eB \ d{w,0) > d{gw,0)}. 

We call l{g) the isometric sphere of g. For a subset S C SL2(C) such that no 
element of S fixes 0, the exterior domain of S is Ext(5') = f]g^gExt(g). The set S 
is a boundary for Ext(S'). A normalized boundary for Ext(S') is a subset S' G S 
such that Ext(S") = Ext(S') and for all g G S' , the geodesic plane l{g) contains a 
face of Ext(S') (i.e. it is a minimal boundary). 

With these definitions it is clear that Do{T) = Ext(F \ {1}). Note that for 
all p G B with trivial stabilizer in F, Dp{r) — uDo{u^^Tu) where u G PSL2(F) 
is such that p = u • 0, so there is no harm in restricting to the Dirichlet domain 
centered at 0. Consider an element g e SL2(C) and A,B,C,D as in formula 
Then 5 • = if and only if C = and, if g does not fix 0, then a simple but lengthy 
computation reveals that l{g) is the intersection of B and the Euclidean sphere of 
center w and radius r, where 

(2) w = ~C-^D and r = 2/|C|, 

and that lnt{g) is the interior of this sphere (the details are in |PaglO[ proposition 
3.1.6]). 

Another property of Dirichlet domains is their rich structure: it gives a presen- 
tation for the group, and also necessary and suflacient conditions for an exterior 
domain to be a fundamental fomain. Suppose F is a Kleinian group in which has 
trivial stabilizer, and let g,h G T. Then we have 1(g) = 1(h) if and only ii g = h. 
We also have gl(g) = I(f;~^), and a point x S 1(g) is in the boundary of Do(T) if 
and only if gx G l(g~^) is too. 

From this, we can group the faces of Do(T) in pairs, one contained in some 1(g) 
and the other contained in l(g~^), and g,g^^ send the faces to each other. This 
is the face pairing structure, and the elements g such that 1(g) contains a face 
of Do(T) are called the face pairing transformations. They generate the group F. 

Now we are going to look for relations. The first type comes from edge cycles: 
consider an edge ei of Do(T) contained in some 1(g) n 1(h), and let gi — g. We 
define inductively a sequence of edges and elements in F in the following way. We 
let e„+i = gnSn- Then Cn+i is contained in 1(5^"'^) n l(gn+i) for a unique I(g„+i) 
(see Figure [TTT|) . If Do(T) has finitely many faces, then the sequence (e„,5„)„ is 
periodic, let m be its period. The sequence of edges C — (ei, . . . , e™) is a cycle of 
egdes, and m is its length. The cycle transformation at ei is ft. = gmgm-i • ■ - ffi, 
and it fixes ei pointwise (property (i)), so it satisfies the cycle relation h'^ = 1 for 
some integer f. 11 v ^ 1, the cycle is called elliptic. At every edge e^, the geodesic 
planes l(g~^) and I(gi+i) make an angle a(ei) inside DqIT). The cycle angle of C 
is a(C) = X^illi ci^(ei)- Since the translates of Dq{T) have to cover a neighborhood 
of ei, we have a(C) = ^ (property (ii)). 

The second type of relations comes from elements of order 2: it may happen 
that 1(g) = l(g~^), then the element g satisfies the reflection relation g^ = 1. 
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Figure 1.1. A length three cycle in a planar cut 



Theorem 3 (Poincare). Let D = Do{T) be the Dirichlet domain of a geometrically 
finite Kleinian group T. Then the face pairing transformation generate the group T, 
and the reflection relations together with the cycle relations form a complete set of 
relations for T. 

Remark 4. In the presentation given by the theorem we consider only one element 
for each pair of face-pairing transformation g,g^^. If we take both in the set of 
generators, we have to add the "inverse" relation g g^^ = 1. 

We are now looking for sufficient conditions for an exterior domain to be a 
fundamental domain. There is another necessary condition, coming from cycles of 
some special points at infinity. A point z e dB is a tangency vertex if it is a point of 
tangency z = fnf of two faces / C l{g)J' C 1(5') of Do{T). If zi = I(go)nl(gi) is a 
tangency vertex, then we define a sequence by letting z^+i — gi-Zi — l{g~^)nl{gi+i) 
while Zi+i is a tangency vertex (otherwise the sequence ends at Zi). If such a 
sequence (zi) is infinite and Dq{T) has finitely many faces, then it is periodic. 
Let m be its period; then (zi, . . . , z„i) is a tangency vertex cycle and the tangency 
vertex transformation is h = gmgm-i ■■■91- The fact that B/T is complete implies 
that h is parabolic (property (iii)). 

Actually all these definitions can make sense for any exterior domain. Sup- 
pose Ext(S') is an exterior domain with S' C F a finite normalized boundary. We 
say that it has a face pairing if S* = 5"^ and for every g £ S the image by g of 
the face contained in 1(17) is the face contained in l{g^^) (equivalently, the image of 
every edge of Ext(5') by the pairing transformation of an adjacent face is an edge 
of Ext(S')). This implies that every cycle is well-defined. We say that it satisfies 
the cycle condition if every cycle satisfies the properties (i) and (ii), and that it 
is complete if every tangency vertex cycle satisfies the property (iii). 

Theorem 5 (Poincare). Let D — Ext(S') he an exterior domain with S finite. 
Suppose D has a face pairing, satisfies the cycle condition, and is complete. Let F' 
be the group generated by the face pairing transformations. Then D is a fundamental 
polyhedron for F'. 

Proof. Both theorems are a special case of the second Theorem in |Mas71j . □ 

1.5. Quaternion algebras and arithmetic Kleinian groups. We can now de- 
scribe the construction of arithmetic Kleinian groups, using orders in quaternion 
algebras. The reader can find more about quaternion algebra in |Vig80| . A quater- 
nion algebra B over a field -F is a central simple algebra of dimension 4 over F. 
Equivalently, if char F ^2, there exists a, 6 G F^ such that B = F + Fi + Fj + Fij 
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with multiphcation given by P ~ a, p — b, ji — —ij; such an algebra is 
written B = (^^^ ■ A quaternion algebra is either isomorphic to the matrix 

ring (-f), or a division algebra. Given an element w ~ x+yi+zj+tij G (^^^ , we 

define its conjugate w = x — yi — zj — tij, its reduced trace trd(w) — w + ili = 2x £ F 
and its reduced norm nid^w) = ww = x^ — ay^ — bz^ + abt^ G F. 

Let be a number field, let be its ring of integers and let B be a quaternion 
algebra over F. An order O C B is a finitely generated Z^-submodule with FO = 
B which is also a subring. We write Of C the subgroup of elements of reduced 
norm 1. 

A place w of -F is split or ramified depending on whether B ®p F^ is isomorphic 
to the matrix ring or not. The set of ramified places is finite and the discriminant 
of B is the product of the ramified finite places, viewed as an ideal m.'Lp. The 
number field F is almost totally real (or ATR) if it has exactly one complex place. 
A quaternion algebra over an ATR field is Kleinian if it is ramified at every real 
place. 

Theorem 6. Let F be an ATR number field of degree n, B a Kleinian quaternion 
algebra over F and O be an order in B. Let p : B ^ A^2(C) be an algebra 
homomorphism extending a complex embedding of F. Then the group T(0) = 
p{Oi )/{+!} C PSL2(C) is a Kleinian group. It has finite covolume, and it is 
cocompact if and only if B is a division algebra. Furthermore, if O is maximal, we 
have 

(3) Oo™,(r(.„ . 

where Ap is the discriminant of F, (^p is the Dedekind zeta function of F, Ab is 
the discriminant of B and = N{m) ■ Hpi^j^ (l - A^(p)"^) for every ideal 91 

ofF. 

Proof This theorem can be found in |MR03[ Theorems 8.2.2, 8.2.3 and 11.1.3]. □ 

An arithmetic Kleinian group is a Kleinian group that is commensurable with a 
group r(C') as in the previous theorem. The object of the next section is to describe 
an algorithm which, given such a group, computes a fundamental domain for r(C'), 
and a presentation with a computable isomorphism. 

2. Algorithms 

We describe every algorithm in ideal arithmetic. In section [2.51 we explain how 
to actually implement these algorithms using floating-point arithmetic. 

2.1. Algorithms for polyhedra in the hyperbolic 3-space. We start with low- 
level algorithms for dealing with hyperbolic polyhedra. A point in B is represented 
by a vector in C -I- Kj; a geodesic plane not containing is represented by the 
Euclidean center and radius of the corresponding Euclidean sphere; a geodesic not 
containing is represented by the Euclidean center and radius of a Euclidean sphere 
and a basis of a Euclidean plane containing the center of the sphere, such that the 
geodesic is the intersection of B, this sphere and this plane. 

Using these representations, it is an exercise in computational geometry to see 
that we can compute the faces, edges and vertices of a convex polyhedron given by 
a finite set of half-spaces containing (the details can be found in |PaglO[ section 
II. 3. 3]). A harder task is to compute the volume of such a polyhedron. We describe 
an algorithm here; it is essentially the same as the one described in |MR03i section 
1.7] but for the sake of completeness we provide all the details here. 
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Algorithm [T] computes the volume of a convex polyhedron with finitely many 
faces. 



Algorithm 1 Volume of a convex polyhedron 

Input: A convex polyhedron P with finitely many faces 
Output: The hyperbolic volume of P 

1: Split every face of P into triangles 

2: Split P into tetrahedra 

3: Using the map rj^^, send every tetrahedron back to H'^ 

4: Express every tetrahedron as a difference of two tetrahedra, each having a 

vertex in the sphere at infinity 
5: For every tetrahedron having a vertex in the sphere at infinity, apply an isom- 

etry to map it to a tetrahedron with one vertex at oo and the other vertices on 

the unit hemisphere 

6: Express every such tetrahedron as a sum and difference of tetrahedra of the 

same type having one vertex at j 
7: Express every such tetrahedron as a sum and difference of tetrahedra of the 

same type with projected Euclidean triangle having a right angle not at 
8: For every such tetrahedron, compute the angles a and 7 and use Proposition [2] 

to compute the volume 
9: Vol(P) -s— sum of every contribution 
10: return Vol(P) 




Figure 2.1. Step g] in Algorithm [T] 

Remarks 7. 

• For step [TJ choose a vertex of the face and link it to every other vertex; 

• For step [21 choose a vertex of P and link it to every computed triangle; 

• For steplH choose an edge and consider a geodesic ray containing it, then 
the tetrahedron appears as the difference between two tetrahedra, each 
having the geodesic ray as an edge and a face of the initial tetrahedron as 
a base (see Figure I^TTt : 

• In step [HI the signs that appear in the sum are the signs of certain deter- 
minants; 
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• In step [51 the angle a is an angle in a Euclidean triangle and can be com- 
puted by elementary trigonometry, and since the upper half-space model is 
conformal, the angle 7 is the Euclidean angle of intersection of the sphere 
and plane representing the faces of the tetrahedron. 

The values of the Lobachevsky function are computed with the following lemma. 
It may be well-known, but we include it for the sake of completeness. 



Lemma 8. For all 9 G (— tt, vr) we have the formula 

61 I 3 - In ^ 



£(61) =7rln 
and the bounds 



C(2n) 



■(- 

^ n(2n + 1) Vtt 



2\e\ 1- 



2 

< - 
- 3 



C(2n) -1/6* 



V 

^ n(2n + 1) Vtt 

n— 1 ^ ^ 



2n ' 



1 



E 



C(2n) 



i{2n + 1) Vtt 



2n 



< 



1 - 



2n 



2tt 



2r+2 



2r+2 



Proof. To derive the first expression we use the previous power series expansion and 
extract the first term of the series expansion of the zeta function. For all 9 G (— tt, tt) 
we have 



^ n(2n+l) \n) ~ ^ n(2n + 1) WJ ^ ^ n(2n-f 1) \n) 

since all these series converge. We only need to compute the power series that 
appears. By derivating twice one finds that for all x G (—1, 1) we have 

°° 2,2n+l /l-X 

V— = 2a;-2:ln(l-a;^)+ln( 

n—l ^ ' 

Letting x = in this expression gives the first formula. 

To prove the inequalities we are going to bound the values C(2'^) and C(2'^) — 1 
for n > 1. By series-integral comparison we get 



< 



2n- 1 

/c— r 

which gives 

C(2n) = H- V fc"^" < 1 + < 2 

for the first value, and 



2n- 1 

fc=2 



C(2n) - 1 = 2-2" + y yt-2n < (i + ')2-2« <3.2 

^-^ \ 2n — 1 ' 



fe=3 



l-2n 



for the second one. Using these inequalities and the bound ^{'2^+1) — ^' 
computing the geometric sum gives the result. □ 
Remarks 9. 

• With the same method, for any k one can obtain a formula with remainder 
termO((^n. 
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• In practise, we precompute the coefficients of the power series we are using. 
By periodicity and oddness, we can always reduce to the case where 9 G 
[0,^]: if the precision is fixed, we know a priori the maximal number of 
terms needed to evaluate the Lobachevsky function. 

2.2. The reduction algorithm. When we have a fundamental domain, it is nat- 
ural to ask for an algorithm which, given any point in the hyperbolic 3-space, 
computes an equivalent point in the fundamental domain and an element in the 
group that sends one to the other. 

Definition 10. Let S* be a subset of a Kleinian group F. A point z E B is S-reduced 
if for ah g e S, we have d(z, 0) < d{gz, 0), i.e. if z G Ext(S'). 



Algorithm 2 Reduction algorithm 

Input: A point w G B, a, finite ordered subset S C PSL2(C) 

Output: A point w' and an clement S E (S) s.t. w' is S'-reduced and w' — 6w 

1: w' ^ w, S ^ 1 

2: g^l 

3: repeat 

4: w' <— gw' , S ^ gS 

5: g ^ the first g E S such that d{gw' , 0) is minimal 
6: until d{gw',0) > d(w',0) 
7: return w',6 



Proposition 11. Given S a finite subset of a Kleinian group T and a point w E B, 
Algorithm\^ returns a point w' and 5 E (S) such that w' is S-reduced and w' — Sw. 

Proof. After step 21 we have w' — Sw and S E {S). Because of the loop condition, 
while the algorithm runs the distance d{w' , 0) decreases. But w' stays in the F-orbit 
of w' and this orbit is discrete, so the algorithm terminates, and when this happens, 
g is an element in S such that d{gw',Q) is minimal and d{gw',0) > d(u;',0), so w' 
is S'-reduced. □ 

Remark 12. At stcpO the g achieving the minimal d{gw' , 0) may not be unique. 
We can then pick any of these elements. Ordering S gives us a canonical choice. 

Reducing points can give interesting information about the elements of the group, 
because if w has a trivial stabilizer, then the orbit map 7 i— > 7 • w is a bijection. 
This is the reason for introducing the following definition: 

Definition 13. Let S" be a subset of a Kleinian group F and w E B. An element 7 E 
PSL2(C) is {S,w)-reduced if 7W is S'-reduced, i.e. if E Ext(S'). 

Given a finite S, w and 7, we can now compute an (S', K7)-reduced element 7 
such that 7 = 7 (mod S) as follows: we reduce with respect to S'; if i5 G {S) 
is such that 5{'^w) is S'-reduced, then 7 = i57 is (S', w) -reduced. We also write the 
reduced element 7 = Reds(7;u') and simply Reds(7) = R.ed5(7;0). A priori this 
reduced element could depend on the chosen ordering in Algorithm [21 

Proposition 14. Suppose that Ext(S') is a fundamental domain for (S). Then 
for w E B outside of a zero measure, closed subset of B, the following holds: for 
every 7 G F, there exists a unique {S,w)- reduced 7 = 7 (mod S). If w E Ext(S) 
then ^ — \ if and only if ^ E (S) . 
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Proof. Let w e r-Ext(5'). The existence follows from Algorithm [51 For uniqueness, 
suppose 7 and 7' are (S", w)-reduced and 7 = 7' = 7 (mod S). Then jw,^'w € 
Ext(5'), and since w is in the orbit of Ext(5'), they are in fact in Ext(5'). But these 
two points are in the same (iS')-orbit, so 7 = 7'. Now assume w S Ext(5'). If 7 = 1 
then 7 = 7=1 (mod S), i.e. 7 G (S). If 7 G (S) then 7=1 (mod S) and 1 is 
(5, ui) -reduced so by uniqueness 7=1. □ 

Since this provides an algorithm to write any element of the group as a word in 
the generators and to compute modulo (S) (with explicit unique representatives), 
this particular kind of generating set deserves a name. 

Definition 15. A subset S* of a Kleinian group F is a basis if Ext(5) is a funda- 
mental domain for (S) — T. If S* is also a normalized boundary for Ext(S'), it is 
called a normalized basis for F. 

2.3. Normalized basis algorithms. Now we describe a general algorithm that 
computes a normalized basis for a cocompact Kleinian group F, which we will 
then apply to arithmetic groups. First note that, after conjugating the group by 
a suitable element in PSL2(C), we may assume that E B has a trivial stabilizer 
in F and that every elliptic cycle has length 1. 

We will use two blackbox subalgorithms. Enumerate and IsFullGroup: 

• Enumerate(F, n) takes as an input a positive integer n and returns a finite 
set of elements in F (the integer n is a parameter for iteration, it does not 
have any mathematical meaning); 

• IsFullGroup(F, S) takes as an input a finite normalized basis S for a sub- 
group (S) C F and returns true or false according to whether (S) = F or 
not. 

In every algorithm, an exterior domain Ext(S') with finite S is represented as a 
polyhedron in B. We begin with a naive algorithm. 



Algorithm 3 Naive normalized basis algorithm 

Input: A Kleinian group F 
Output: A normalized basis S for F 

1: S ^9, n^O 

2: repeat 

3: repeat 

4: n n + 1 

5: add Enumerate (F, n) to S 

6: S* normalized boundary of Ext(S') 

7: until Ext(5') has a face-pairing and Ext(S') is complete and Ext(5) satisfies 

the cycle condition 
8: until IsFuUGroup(F, S) 
9: return S 



We say that Enumerate is a complete enumeration of F if we have 



Proposition 16. // F is geometrically finite and Enumerate is a complete enu- 
meration of F, then Algorithm\^ terminates after a finite number of steps and the 
output S is a normalized basis for F. 




n>0 
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Proof. The Dirichlet domain centered at for F has finitely many faces by geo- 
metric finiteness. Since Enumerate is a complete enumeration, a boundary for this 
Dirichlet domain will be enumerated after a finite number of steps. The algorithm 
will then terminate as all the conditions are satisfied by Dirichlet domains. The 
output will then be a normalized basis for F by Step |B] and Theorem [S] □ 

We will now use the reduction algorithm to improve upon Algorithm [3l The 
main ideas are 

• reducing the elements that we have to find smaller ones 

• when the face-pairing condition, the cycle condition or the completeness 
condition fails, using this fact to find elements that make the exterior do- 
main smaller. 

For clarity, we divide Algorithm S] into four routines. Algorithm 0] uses these rou- 
tines to compute a normalized basis for a geometrically finite Kleinian group F. 



Algorithm 4 Normalized basis algorithm 

Input: A Kleinian group F 

Output: A normalized basis S for F 

1: ^ 0, n^O 

2: repeat 

3: repeat 

4: n n -t- 1 

5: add Enumerate (F, n) to S 

6: S <- KeepSameGroup(S') 

7: S ^ CheckPairing(S') 

8: 5 ^ CheckCyclcCondition(S') 

9: S ^ CheckComplete(S') 

10: until Ext(S') does not change 

11: until IsFullGroup(F, S) 

12: return S 



The first routine, KeepSameGroup, reduces elements as much as possible to 
eliminate redundant ones and find smaller ones. 



Algorithm 5 KeepSameGroup 
Input: A finite subset S C PSL2(C) 

Output: A new S generating the same group with smaller elements 
1: repeat 

2: [/ normalized boundary of Ext(S') 

3: for all g e S do 

4: g-i^Reduig) 

5: if 5 7^ ±1 then 

6: add g to U 

7: end if 

8: end for 

9: S 

10: until Ext (5) does not change 
11: return S 



Proposition 17. If S is a subset of a Kleinian group, then Algorithm\^ terminates 
and does not change the group generated by S. 
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Proof. We first prove the second claim. Every element added to S belongs to the 
group generated by S as it is a reduction by J7 C S" of an element in S. Moreover, 
every element that is discarded has Red;7(5) — ±1 so at the end of the loop we 
have g G {S), and every other element g £ S\U is replaced hy g = Redu{g) & {U)g, 
so the group generated by S does not change. 

Now we prove that the algorithm terminates. First consider the initial S. 
Let M = max{d(5 • 0,0) : g € 5} and Xq = {g e (S) : d{g -0,0) < M}. The 
set Xq is finite since (S) is a Kleinian group, and we have S C Xq. By definition of 
reduction, every element added to U is in Xq. Moreover, by Step [2] if an element g 
is discarded then its isometric sphere l{g) does not intersect Ext(S'), so g • is in 
the complement of Ext(S'): g cannot be the reduction of any element, so it cannot 
be added again. Similarly if g e 5 \ C/ is replaced hy g g, then g is not reduced 
so it cannot be added again. Hence the algorithm terminates. □ 

The second routine, CheckPairing, checks whether Ext(iS') has a face-pairing. If 
it does not, it finds elements that make Ext(S') smaller. 



Algorithm 6 CheckPairing 
Input: A finite subset S C PSL2(C) 

Output: A new S such that Ext(S') is smaller if it did not have a face-pairing 
1: S^S\JS-^ 

2: for all e edge in l{g) and g G S, s.t. ge not an edge of Ext(S') do 

3: X -(^ X & e such that gx ^ Ext(S') 

4: g^Reds{g;x) 

5: add g, g~^ to S 

6: end for 

7: return S 



Proposition 18. //Ext(S') does not have a face-pairing, then after applying Algo- 
rithm\^ Ext(5') is strictly smaller. 

Proof. If there is a nonpaired edge, at Step [5l since x £ l{g) we have d{gx, 0) = 
d(x, 0) and since gx ^ Ext(S') we have d{gx,0) > d{gx,0). Putting these two 
together gives d{gx,0) < d(x, 0), i.e. x € Int(5) so finally we have Ext(S'U {g}) ^ 
ExtCS*). ' □ 

We give a second possible algorithm for CheckPairing, which is simpler but 
less efficient in practice. It uses the fact that if a non-elliptic cycle has length 
three (which is generically the case), then it is of the form e C 1(g) n l{h), ge C 
I(ff-i) n ligh-^), he C l{hg-^) n I(/i-i). 



Algorithm 7 CheckPairing' 
Input: A finite subset S C PSL2(C) 

Output: A new S such that Ext(S') is smaller if it did not have a face-pairing 
1: S^SUS-^ 

2: for all g,he S s.t. l{g) (1 l{h) ^ 9 and h ^ g-^ do 
3: add gh^^,hg^^ to S 
4: end for 
5: return S 



Proposition 19. //Ext(S') does not have a face-pairing, then after applying Algo- 
rithm^ Ext(5) is strictly smaller. 
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Proof. If there is a nonpaired edge, then there exists elements g,h G S in the 
normahzed boundary of Ext(S') and a point x e 1(5^^) n Ext(S') such that g^^x S 
Int(/i) (so that h ^ g'^). Since we also have g ^x G 1(g) and l{g) is not contained 
in Int(/i), we get 1(5) n l{h) ^ 0. But then d(a;, 0) = d{g'''^x,0) > d{g~'^hx,0), 
so xe lnt{g-^h): we have Ext(S' U {g^^h}) £ Ext(S'). □ 

Remark 20. Although this algorithm is less efficient than Algorithm [SJ it is in- 
teresting as it gives a geometric understanding of the method described in |Lip02| : 
"we consider words that are two-word combinations of those forming the sides of 
the existing domain to modify the domain. (...) This procedure has proven to be 
fast and effective in practice." Proposition 1191 explains why taking products of two 
elements forming the sides of the domain is useful, and in Algorithm [7] we get a 
geometric description of the the products that one should form. Actually, the com- 
putation in the proof of Proposition [7] also shows that if l{g~^h) reduces Ext(S'), 
then 1(g) n l{h) ^ 0. 

The third routine, CheckCycleCondition, checks whether Ext(S') satisfies the 
cycle condition. If it does not, it finds elements that make Ext(5) smaller. 

Algorithm 8 CheckCycleCondition 
Input: A finite subset S C PSL2(C) 

Output: A new S s.t. Ext(S') is smaller if it did not satisfy the cycle condition 
1: Compute every well-defined edge cycle 
2: for all g cycle transformation for the edge e do 
3: if g 7^ ±1 fixes at most one point in e then 
4; 5^5U{5,5-H 

5: else if g 7^ ±1 fixes every point in e then 

6: s^sy^ (g) 

7: else 

8: m length of the cycle 

9: for all < i < TO do 
10: gi...gi 
11: add h, h^^ to S 

12: end for 
13: end if 
14: end for 
15: return S 

Remarks 21. 

• If we assume that every non-elliptic cycle has length three, then the steps |H1- 
[T^ are unnecessary, as in this case the partial cycle transformations at an 
edge contained in l{g) n l{h) are g, h = {hg~^)g, 1 = h~^{hg~^)g. 

• If we know in advance that the group F is torsion-free, then we can omit 
the steps [5HS1 

• Assuming both, we can omit CheckCycleCondition completely. 

Lemma 22. Suppose S d T is a subset of a Kleinian group F such that has 
a trivial stabilizer in T, and suppose there is an element G F \ {±1} and a 
point X e Ext(S') such that hx € Ext(S'). Then Ext(5 U {h, h^^}) $1 Ext(S'). 

Proof. First suppose that d(a::, 0) < d{hx,0). Then writing x — h^^{hx) — h~^y 
we get d{h~^y,0) < d{y,0) i.e. y e Int(/i). But we also have y e Ext(S') so Ext(S'U 
{h}) g Ext(S'). 

Othewise we have d{hx,0) < d(a;,0). This means that x G lnt{h^^), but 
since x G Ext(S') we get Ext(S' U {/i^H) £ Ext(5). □ 
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Proposition 23. //Ext(S') does not satisfy the cycle condition, then after applying 
Algorithm\^ Ext (5) is strictly smaller. 

Proof. Since the cycle transformation at an edge stabilizes it, if the edge is not equal 
to a geodesic then the cycle transformation fixes it pointwise and condition (i) is 
automatically satisfied. Suppose that there is a cycle for an edge e equal to a 
geodesic and that does not satisfy condition (i), and let g be the corresponding 
cycle transformation. Then the transformation g is either loxodromic, or elliptic of 
order 2 with exactly one fixed point in e. In both cases, Step 2] is executed. In the 
first case, since the interior of the isometric sphere of a loxodromic element contains 
one of its fixed points and the interior of the isometric sphere of its inverse contains 
the other, we have Ext({g,g"^}) n e £ e so Ext(S' U {g,g~^}) £ Ext(S'). In the 
second case, the edge e contains exactly one fixed point of g in H^, so we again 
have Ext({5}) n e £ e and we get Ext{S U {g, g'^}) £ Ext(S'). 

Now suppose some cycle angle for a non-elliptic cycle is larger than 27r. Then 
considering the images P = Ext(5), . . . ,{gi . . .gi)^^P oi P — Ext(5) that 

glue one after another around e, there is an overlap: there exists a point x G P 
such that hx £ P for some h considered in Step [TU] But then after Step [TT] we 
have Ext(S' U {/i, h^^}) Ext(5) by Lemma [22] Since the cycle transformation is 
the identity, the angle cannot be smaller than 27r. 

Finally suppose some cycle angle for an elliptic cycle at an edge e with cycle 
transformation g with order v does not satisfy condition (ii). The cycle has length 1, 
so e C l{g) n 1(5^^), and the angle at e is a multiple of But after Step [S] 
Ext({(7, (7~^}) is replaced by the Dirichlet domain of the finite group (5), which 
satisfies the cycle condition, so the new angle at e is equal to □ 

The fourth routine, CheckComplete, checks whether Ext(S') is complete. If it is 
not, it finds elements that make Ext(S') smaller. 



Algorithm 9 CheckComplete 
Input: A finite subset S C PSL2(C) 

Output: A new S such that Ext(S') is smaller if it was not complete 
1: Compute every tangency vertex cycle 
2: for all h tangency vertex transformation do 
3: if /i 7^ 1 is loxodromic then 
4: add h, h^^ to S 
5: end if 
6: end for 
7: return S 



Remark 24. If we know in advance that the group F is cocompact, we can omit 
CheckComplete in Algorithm 0] and simply test whether Ext(5) is bounded. 

Proposition 25. // Ext(S') is not complete, then after applying Algorithm 
Ext(S') is strictly smaller. 

Proof. If /i is a tangency vertex transformation at 2; = 1(g) n l{g') G dB, then 
it fixes z. By looking at the successive images of the polyhedron along the cycle 
one sees that l{g') separates l{g) from hl{g), so h has infinite order. If Ext(S') 
is not complete, then h is loxodromic. But then z E Int(/i) U lnt{h~^) so we 
get Ext(5 U {/i, /i-i}) g Ext(S'). □ 

Proposition 26. Let T be a Kleinian group. The following holds for Algorithm^ 
applied to T: 
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(i) Suppose the algorithm terminates. Then the output is a normalized basis for T. 
(ii) Suppose that V is geometrically finite and Enumerate is a complete enumera- 
tion ofT. Then the algorithm terminates. 

Remark 27. In practise Algorithm U] runs much faster that the naive Algorithm [3] 
(see section 13. l.ip . but unfortunately we could not prove it. What we believe is that 
in Algorithm U] the blackbox Enumerate only needs to find a set of generators for 
the group, and then the other routines find the elements of the normalized basis; 
in Algorithm [3] the blackbox Enumerate needs to find directly the elements of the 
normalized basis, which is harder. The natural idea would be to put the routines 
in a loop that would not contain Enumerate in Algorithm 21 but then it is not clear 
whether this internal loop would terminate; actually in general it is false, since T 
may admit finitely generated subgroups that are not geometrically finite. 

Proof. 

(i) If the algorithm terminates, then by Theorem[Sl since Ext(S') is complete, has 
a face-pairing and satisfies the cycle condition, the set 5 is a normalized basis 
for {S). It is then valid to use IsFuUGroup to check that {S) = T. 

(ii) The Dirichlet domain centered at for F has finitely many faces by geometric 
finiteness. Since Enumerate is a complete enumeration, a boundary for this 
Dirichlet domain will be enumerated after a finite number of steps. The 
algorithm will then terminate as all the conditions are satisfied by the Dirichlet 
domain. 

□ 

2.4. Instantiation of the blackboxes. 

2.4.1. Enumerate and IsFullgroup for a group given by generators. If the group T is 
given by a finite set of generators G. We can take for Enumerate the algorithm that 
writes every word of length n in the generators, and we can take for IsFuUGroup 
the algorithm that reduces every element in G with respect to the given normalized 
basis S and returns whether every generator reduces to ±1: by Proposition [TJl this 
is equivalent to F C {S). 

2.4.2. Enumerate and IsFullgroup for an arithmetic group. We provide a possi- 
ble instantiation of the blackboxes Enumerate and IsFullgroup for an arithmetic 
group F(C') attached to a maximal order O in a Klcinian quaternion algebra B 
with base field F of degree n. 

We describe IsFullgroup first. A subgroup is proper if and only if it has a 
(possibly infinite) covolume of at least twice the covolume of F (the quotient of the 
covolumes is the index of the subgroup). Since F comes from a maximal order, 
the covolume of F is given by (jS]), which we can compute, and the covolume of a 
subgroup can be computed with Algorithm [T] once we have a normalized basis. We 
take for IsFullgroup the algorithm that computes the covolume Covol(F) by the 
formula and the volume V of Ext (5) for the given normalized basis S, and returns 
whether ^-.^^^^^^ < 2. Since S' is a normalized basis for {S), the polyhedron Ext (5) 
is a fundamental domain for {S) so the volume V equals the covolume of {S). 

We now describe an instantiation of the blackbox Enumerate for the Kleinian 
group associated with an order O in _B. Under the natural embedding O C B ^ 
B (S)Q K, the order O is discrete. Now suppose that we have a positive definite 
quadratic form Q : B M.. Then O becomes a full lattice in a real vector space 
of dimension An. We can use lattice enumeration algorithms such as the Kannan- 
Fincke-Pohst algorithm jFP851 IKan83| to enumerate elements in O that are short 
with respect to Q. We can then select the elements having reduced norm 1. As we 
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increase the bound on the values of Q, we will get every element in Oi . A priori 
any such quadratic form would work, but here we describe one that has a geometric 
meaning. 

Recall we can embed B in M2 (C) so that Of becomes discrete in SL2 (C) . This 
embedding is only defined up to conjugation by an element of PSL2(C). Let p be 

such an embedding. If S = (^^^ we can take for example 

/ X + ya z + ta 
P--'' + y' + '-^+^'^^[iz-ta)P x-ya 

where ct is a complex embedding of F, /3 — a{b) and a is a square root of a{a). 
For m— (° ^) G A^2(C), we define invrad(r7i) = |(c + 5) + (d— a)j\ . 

Proposition 28. The quadratic form Q : S (g) K. — > M defined for all x E B by 

Q{x) = invrad(p(x)) + tri;-/Q(nrd(a;)) 

is positive definite and satisfies 

Qi-)-^^^^+-forallxeO^ 

where rad(g) denotes the Euclidean radius of the isometric sphere of g £ SL2(C) 
3 ■ 7^ Oi '^'^'^ °o otherwise. 

Proof. We show first that Q is positive definite. For a matrix m G A^2(C) we 
have invrad(m) = |c + 6p + \d — = ||m|p — 25R(detrn) where || • || is the usual 
norm on A^2(C), so that || • |p is a positive definite quadratic form on A^2(C). 
Since nrd is a positive definite quadratic form on H and we have the decomposi- 
tion B (g) R = Al2(C) ® H""^, we can construct a positive definite quadratic form 
on i? (Ki M by letting for all x e B ®R 

Q{x) = ||m||^ + nrd(/ii) + h nrd(/i„_2) = invrad(TO) + tr^^R/R(nrd(x)) 

where 

x^m + hi + --- + h.n-2 e M2{C) ® H"-2^ 

since 2 5R(detm) + nrd(ft.i) + • • • + nrd(/i„_2) = tr^^R/R(nrd(x)). This gives the 
positive definiteness. 

For the formula on Of , note that according to it is 

invrad(g) = | (c + 6) + (d - a)i\ = ^^^^^^^^ 
for g G SL2(C) not fixing in S, and if g fixes then invrad(5) =0. □ 



We obtain the following enumeration algorithm, which is a complete enumeration 
of r(0). It depends on a parameter: a sequence of bounds An — > 00. 



Algorithm 10 Enumerate 

Input: A positive integer n 
Output: A finite subset L C r(0) 
1: L ^ 

2: for all x e O such that Q{x) < An do 

3: if nrd(a;) = 1 then 

4: Add p{x) to L 

5: end if 

6: end for 

7: return L 



17 



We are now going to present a non-deterministic enumeration algorithm. It is 
not a complete enumeration, but performs better in pratice (see section r3.1.2p . It 
uses variants of the former quadratic form. 

Definition 29. Let zi,Z2 G T-l^- Let /ii,/i2 G SL2(C) be such that zi = hi ■ j 
and Z2 = /i2 • j- We then define the quadratic form Qzi.z2 by 

Qzi,z2ix) = invrad(/i2"V(2;)/ii) + tri^/Q(nrd(a;)) 

for all X € B. 

This family of quadratic forms has the following properties. 

Proposition 30. Let zi,Z2 G ■ Then Qzi.zz does not depend on the choice 
of /ii, /i2 G SL2(C) such that zi ^ hi ■ j and ^2 = ^2 ■ J- It is positive definite, and 
for all g £ Of we have 

Qzuz2{9) = 2coshd(5Zi,Z2) - 2 + n. 

Proof. The matrices hi and /12 are defined up to right multiplication by SU2(C) 
(the stabilizer of the point j). For all matrices m e A^2(C) we have invrad(TO) = 
||to||^ — 2 3fi(detm), which is not changed by left and right multiplication of m by 
elements of SU2(C), so that Qzi,z2 does not depend on the choice of hi and ft.2- 

For zi = Z2 — j the formula reads \\g\\^ — 2coshd{gj,j) for all g G SL2(C), 
which is well-known (and is a direct consequence of the explicit formulas for the 
hyperbolic distance). But then for arbitrary zi, Z2 G we have 

= 2 cosh d{h2 ^g hi = 2coshd{ghij,h2j) = 2 coshd(5Zi, 22)- 

□ 

This family of quadratic forms is very useful, as it enables us to determine the 
elements g G r(C') such that gzi is close to Z2. We propose the following non- 
deterministic algorithm for enumerating elements in T(0). It depends on a choice 
of some parameters: an increasing sequence of radii i?„ — >■ 00 (the radius of the 
search space), a sequence of positive integers Nn G Z>o (the number of enumeration 
in small balls) and a bound A (bound on the quadratic form used). For tui, W2 G B, 
we write Qwi,w2 '■— Q-n-'^{wi),rr^(w2)- 



Algoritiim 11 Enumerate' 

Input: An positive integer n 
Output: A finite subset L C T{0) 
1: L ^ 

2: for i = 1 to A„ do 

3: Draw a point w € B such that d(0,u') < Rn randomly, uniformly w.r.t. the 

hyperbolic volume 
4: for all X E O such that Qo.w{x) < A do 
5: if nrd(x) = 1 then 
6: Add p{x) to L 

7: end if 
8: end for 
9: end for 
10: return L 

Remarks 31. 

• We can also use these quadratic forms differently: if we miss an element of 
the group to "close off" the exterior domain around a point at infinity ^, 
we can look for elements of small Qj^z where z — ^. This is a similar idea 
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as in Remark 4.9 in |Voi09j . but the quadratic form that was used there is 
the analogue of Qz,z- If 9 is the element that we are looking for, d{gz,z) 
is bounded by below by a positive constant if g is loxodromic (which is the 
generic case), while d(f/j, z) — >• as z — gj. 
• The efficiency of this algorithm depends on the choice of the parameters 
R„ and A. Heuristics led us to the following choice, which works well in 
practise: 

— we use a small bound A — a ■ \ ApN{Ab)\ "i-^-'™ so that the number 
oix € O such that Qa,w{x) < yl is approximately constant by Gaussian 
heuristic; 

— experimental evidence and [BGLSlOl Theorem 1.5] suggest that a num- 
ber of random elements of F proportional to Covol(r) has a good prob- 
ability to generate F, and by Gaussian heuristic one needs 0(Covol(F)) 
random centers to obtain one element of the group on average, so we 
choose Nq = (3 ■ Covol(F)^, and we increase it exponentially fast: Nn = 
(l + r/fA^o; 

— the radius i?„ has to be large enough to ensure good randomness of the 
elements of F, so we choose Rq such that Vol(B(w, i?o)) — Covol(F)''' 
and we increase it in arithmetic progression (so the volume increases 
exponentially fast): Rn — Ro + ^ ■ n. Because of our choice of Nn we 
take 7 > 2. 

Now we explain how we draw points at random in the ball B(0, R) of radius R. 
Since the hyperbolic volume is invariant by rotation around 0, it is equivalent to 
draw a random point uniformly on the sphere, and then multiply it by an appro- 
priate random scalar independent from the point on the sphere. Thus we only have 
to determine the distribution of the distance from of the points in the ball of 
radius R. Let X be a random variable with uniform distribution in B(0,i?). The 
cumulative distribution function of the distance to is 

Recall that we have v{r) :— Vol(B(0, r)) — 7r(sinh(2r) — 2r). It is then clear that 
the function fji : [0,i?] — ^ [0, 1] is a continuous bijection. It implies that d(0,X) = 
f^^{U) where U is a uniform random variable in [0, 1]. We rewrite that expression 
as d{0,X) — v^^{U') where U' is a uniform variable in [0,v{R)]. It is well-known 
how to draw a uniform variable in an interval and on a sphere, and can be 
computed by Newton iteration. 

2.5. Floating-point implementation. Here we describe a floating-point imple- 
mentation of the above algorithms. We start with a lemma giving us control on the 
error made when having an element of the group act on a point. We only study 
the stability of the algorithm, so we do not take into account the error made by 
rounding in elementary operations. 

Lemma 32. Let g e SL2(C), g G A^2(C) and w,w G B. Let e — \w—w\, rj = 5I 
and 5 = ■ Suppose that (||g||e + 2?])^ < Then the quantity gw obtained by 

applying Formula ^ to g and w is well-defined, and we have 

\g-w~gw\< 68 S^\g\\^e + 136 Si\\g\\^T]. 

Proof. By direct computation we have |A — < ^/2ri and \A\ < \/2\\g\\, and the 
same inequalities for B,C, D. We write 

g-w = {Aw + B){Cw + £»)-! = ,^ ^ r^^o i^w H- B){wC + D) 

\Lw -\- iJp 
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and similarly for g^w. Another direct computation gives 

4 

\Cw + D\ 



(4) -\9-w\^=[l- ^^-^ ) (hp - 1) 



which shows that 



<-(l + 2^)<-^ and ^^I'^^ + ^l' 



\Cw + D\^ - 4' ' - A 3(5 
By the triangle inequality, adding and substracting Aw gives 

\Aw - Aw\ < V2\\g\\e + V2r] 
and the same inequality for Cw. We get 

\{Cw + D)-{Cw + D)\^ <2{\\g\\e + 27j)^ < '^""^^^'^ 

since by hypothesis we have (||5||e + 2ri)^ < jg. In particular Cw + D and gw 
is well-defined. By the mean value theorem this gives 

\\Cw + D\~^ -\Cw + D\-^\ < (6<5)^(|l.g|le + 277) 

We also get 

\{Aw + B){wC + D) ~ {Aw + B){wC + D)\ 

< \Aw + B\{V2\\g\\e + 2V2f]) + 2\Cw + i:i|(V2||(7||e + 2V2r]) 

< i2V2\\g\\)iV2\\g\\e + 2V2v) + {2V2\\g\\)i2V2\\g\\e + 4V2v) 

12\\g\\\ + 24\\g\\rj. 

Finally we have 

\{Aw + B){Cw + D)-^ ~ {Aw ^- B){Cw + bY^\ 

< \g ■ w\\Cw + D\H66)H\\g\\e + 2r^) + ^^^{I2\\g\\\ + 24||gh) 

< (24^6 + 9)(5i||5||3e + (48^6 + 18)<5i||5||277 

< 68(5i||5|h\ + 136<5i||g|p?7 

as claimed. □ 

In the following, we want to always have (|l(7||e + 2rjY < for every element g 
and every point w considered, where e is the imprecision on the points in B, rj the 
imprecision on the elements g considered, and 77 = |e. 

We now describe the modification of the algorithms for the fioating-point version. 
In the reduction algorithm (Algorithmic]), we choose a > and in Step H] we replace 
the inequality d{gw',0) > d(w;',0) by \ (j^f^jy\-i < 1 + a. Since we have w' G Int((7) 
if and only if \Cw' + -Dp < 4, the modified condition is indeed an approximation 
of the exact condition. 

Proposition 33. Let l3 ^ a - 68 SiAPe ~ 136 Siu'^rj where 5 = yrr^ and M = 
maxggs \\g\\- If (3 > 0, then the floating-point version of the reduction algorithm 
terminates. 

Proof. Formula Q may be rewritten 

which gives, if the modified condition of Step |6] is not satisfied 

l-\g-w'\^>{l + a){l-\w'\'). 

Lemma [31] gives 
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so that 1 — is multiplied by 1 + /3 at each step of the algorithm. Since we also 
have 1 — jw'p < 1, the algorithm terminates. □ 



We want to use a uniform a which tends to as e — >■ 0. For this, we assume 
that we only consider points w such that 1 — |wp > 2e3 and elements g such 
that llgll < e~9. Assuming that e < 10~^ we can then take a = ISes. These 
assumptions also ensure that (Uglje + 277)^ < and are compatible since the 
points g ■ that we have to consider satisfy 1 — |g • Op = ||^||'^_^^ > > 2e^ . 

There is no change in KeepSameGroup (Algorithm [S]) : the same argument shows 
that the algorithm terminates, regardless of finite precision in the computations. 

The routine CheckPairing (Algorithm ^ should only consider an edge e con- 
tained in l{g) as not being paired if there is x g e and we have the stronger 
inequality \Cgx + Dp < and C,D correspond to h for some h & S. This 
ensures that the floating-point reduction will yield a non-trivial element, since at 
least one step of reduction will be performed. 

The routines CheckCycleCondition (Algorithm [S]) and CheckComplete (Algo- 
rithm [S]) contain only finite loops regardless of the use of finite precision, so there 
is no change in them. 

Proposition 34. The floating-point version of the Normalized basis algorithm (Al- 
gorithm^ terminates. 

Proof. By the arguments above, each of the routines terminates. Moreover, because 
of precision restriction we impose ||.g|| < e^s for every element g of the group 
considered in the algorithm, so that only finitely many g can be used. So the 
algorithm terminates. □ 

Of course if the precision chosen is insufficient, the algorithm may terminate 
with an error or a wrong answer, but with Riley's methods |Ril83) . we we can use 
Poincare's theorem with the approximate fundamental domain to prove that the 
computed presentation is correct. Alternatively, we could check the fundamental 
domain algebraically, but this is likely to be time-consuming. 
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2.6. Master algorithm. As a summary, this is our master algorithm for comput- 
ing an arithmetic Kleinian group associated with a maximal order. 



Algorithm 12 Master algorithm 

Input: A maximal order O in a Kleinian quaternion algebra B 
Output: A finitely presented group G, and two computable group homomor- 
phism <j) : G —> r{0) and tjj : r(C') — >■ G, inverse of each other 
1: Choose an embedding p : B ^ A^2(C) s.t. the point has trivial stabilizer in 

the group r(0) = p(On/{±l} 
2: V ^ Covol(r(e')) computed with Formula Q 
3: function IsFullGroup(S') do 
4: compute V = Vol(Ext(S')) with Algorithm [1] 
5: return V < 2V 
6: end function 

7: Enumerate ^ Algorithm [10] or Algorithm [Til 

8: S ■(r- output of the Normalized Basis Algorithm [3] 

9: R ■(r- inverse, cycle and reflection relations from Theorem [3] 

10: G ^ {S\R) 

11: Let (j> : G ^ r(C') be the map that evaluates words in the generators 

12: Let : T{0) — >■ G be the map that writes elements as words in the generators 

using Algorithm [5] 
13: return G, (f), ip 



Remarks 35. 

• If we want to compute the group that is the image of a smaller order, or 
more generally a finite index subgroup F' of the group F given by a maximal 
order, we can compute first a normalized basis for the larger group F, and 
then compute the index by standard coset enumeration techniques. This 
gives the covolume of the smaller group, and even a set of generators for it, 
so we can then apply the same algorithm that is described here. 

• We may also want to compute a maximal group F" in the commensurability 
class of F as described in |Bor81) . Since it is the image in PSL2(C) of 
the normalizer of O in B, we may use the same enumeration techniques. 
The index is given in terms of a class group and a finite quotient of units 
inZp (see [EorSl] for details), which can be computed, so again we get the 
covolume of this larger group, and can apply the same technique. 

3. Examples 

The author has implemented the algorithm described in the previous section in 
the computer system Magma (BCP97) . Our package KleinianGroups is available at 
http : //www . normalesup . org/~page/Recherche/Logiciels/logiciels-en . html 

Here we show some examples of the output of this code. In sections [3.11 and the 
computations are performed on a 1.73 GHz Intel 17 processor with Magma v2.18-4. 
The more extensive computations of sections 13.31 and are run on a 2.5 GHz Intel 
Xeon E5420 processor with Magma v2. 17-12. 

3.1. Comparison between subalgorithms. 

3.1.1. Comparison between the normalized basis algorithms. Consider the ATR sex- 
tic field F of discriminant —92779 generated by an element t such that t^ — t^ — 
2t4-h3t3-i2_2< + l = 0, andlet Zp be its ring of integers. Let B = (^^) be the 
quaternion algebra ramified only at the real places of F. Let O be a maximal order 
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in B (they are all conjugate). The Kleinian group T{0) has covolume 0.3007. . .. 
We compare our algorithm with the naive Algorithm [31 Both need a precomputa- 
tion of 3 seconds for the computation of the coefficients of the Lobachevsky power 
series and 4 seconds for the evaluation of the Dedekind zeta function at 2. Our 
algorithm then computes a Dirichlet domain in 2 seconds, and enumerates 37 ele- 
ments of O, yielding 21 elements of r(C'). The naive algorithm (actually we only 
removed the routine CheckPairing) computes the same Dirichlet domain in 48 sec- 
onds and has to enumerate 16 246 elements of O, yielding 1713 elements of r(0). 
The fundamental domain (Figure l3TT|) has 18 faces and 42 edges. 




Figure 3.1. Dirichlet domain of a Kleinian group over a sextic field 



3.1.2. Comparison between the enumeration algorithms. Consider the ATR number 
field F of degree 8 and discriminant —407793664, generated by an element t such 
that - 4£' + At^ + 2t^ - St" + At^ + 51"^ -2t-l ^ 0, and let be its ring of integers. 
Let B = { ~^F~^ ) b^ quaternion algebra ramified only at the real places of F. 
Let O be a maximal order in B (they are all conjugate). The Kleinian group r(0) 
has covolume 56.509 .... We compare the performance of our algorithm when 
using the enumeration algorithms [10] or [11] With the deterministic enumeration 
algorithm llOl our code computes a fundamental domain in 12 hours and 45 minutes 
(45943 seconds, most of which is enumeration), and enumerates 84 159 799 vectors, 
yielding 1600 group elements. With the probabilistic enumeration algorithm [TT] 
our code computes the same Dirichlet domain in 71 seconds, and only needs to 
enumerate 3511 vectors, yielding 164 group elements. It spends 2 seconds for com- 
puting the value of the zeta function, 16 seconds for enumeration, 3 seconds for the 
routine KeepSameGroup, 40 for CheckPairing and 10 for computing the volume of 
the polyhedron. The fundamental domain (Figure has 202 faces and 582 edges. 

3.2. Relation to previous work. In this section we show how to recover exam- 
ples covered by earlier work with our algorithm. When available, we provide a 
comparison of running times between public implementations and our code. One 
should keep in mind that these are only comparisons of implementations since the 
complexity of the algorithms is usually unknown. 
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Figure 3.2. Dirichlct domain of a Kleinian group over an octic field 



3.2.1. Bianchi groups. Let F be an imaginary quadratic field with ring of inte- 
gers Ijf- Consider the quaternion algebra B — M.2{F) and the maximal or- 
der O = A^2(2f)- Then the group T{0) — PSL2(Zi?) is called a Bianchi group. 
There exists already several programs computing fundamental domains for these 
groups jRahlOl lYaslOj but they only work for Bianchi groups while ours deals 
with general arithmetic Kleinian groups. Table 13.11 gives the running time (in 
seconds) of our Magma package and other public implementations. The first 
three columns correspond to the discriminant of the field, its class number and 
the covolume of PSL2(Zi?). The last four columns display runnig times in sec- 
onds : Bianchi . gp (RahlO) in GP [The 11) which implements Swan's algorithm 
for PSL2(Zi?), our code KleinianGroups computing PSL2(Zi?), the code pro- 
vided by Magma implementing the algorithm of |YaslO| using Voronoi' theory 
for PGL2(Zi^'), and our code for PGL2(Zf). Note that it is not surprising that 
computing PGL2(Zi?) is faster : the group is larger by an index 2, so the covolume 
is twice smaller and our computation is 4 times shorter (see also section [211]). 

3.2.2. Arithmetic Fuchsian groups. Let be a totally real field and B a quaternion 
algebra ramified at every infinite place but one. Let O be an order in B. Then the 
group T{0) = Of /{±1} embeds into PSL2(M), in which it is discrete with finite 
covolume: it is an arithmetic Fuchsian group. Using the action of PSL2(M) on 
the upper half-plane J. Voight |Voi09) was able to compute fundamental domains 
for these groups. Since we have PSL2(K) C PSL2(C), a Fuchsian group can be 
seen as a Kleinian group leaving a geodesic plane stable. Using this we can also 
compute arithmetic Fuchsian groups with our code. Our probabilistic enumeration 
Algorithm [TT] leads to an improvement in high degree. As an example, consider the 
totally real field F with discriminant 9685993193, generated by an element t such 



that t^ - 2t^ - 7t^ + lli*^ + 15i5 - Ibt"^ - lOt^ + 7t^ + 2t-l^0. Let B = (^^j 

with a = -3t^ + 2f + - - 93t4 + 90<2 + 2f - 26 and fc = -1, which is 
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Magma 


KG, PGL2 


-3 


1 


0.169 


0.015 


0.93 


0.43 


0.83 


-15 


2 


3.139 


0.152 


0.92 


0.8 


2.32 


-23 


3 


6.449 


0.176 


1.22 


1.11 


2.06 


-39 


4 


13.80 


2.37 


9.44 


3.05 


4.36 


-47 


5 


19.43 


3.83 


19.9 


5.33 


6.96 


-71 


7 


37.53 


21.6 


36.6 


17.8 


13.2 


-87 


6 


44.72 


25.7 


45.1 


17.3 


16.4 


-95 


8 


57.06 


41.4 


43.8 


33.9 


19.3 


-119 


10 


82.93 


7080. 


137. 


99.5 


25.6 


-167 


11 


132.3 


1545. 


391. 


188. 


80.9 


-199 


9 


148.5 


3840. 


393. 


224. 


92.7 



Table 3.1. Running times for Bianehi groups 



ramified at every real place but one. Let O be a maximal order in B. The Fuch- 
sian group r(0) has coarea 103.67 . . . ; our code computes a fundamental domain 
for this group in 13 minutes (735 seconds). The code provided by Magma and 
implementing the algorithm of |Voi09| computes a fundamental domain for r(0) 
in 1 hour and 10 minutes (4204 seconds). 

3.2.3. The Hamiltonians over 'L \^^i^ -\ . Consider the field F = Q(-\/— 7) and the 
quaternion division algebra . Then O = 'Lp + Zj^i + Zj^j + Z^^ij is a 

non-maximal order in _B. A fundamental domain for this group was computed by 
C. Corrales, E. Jespers, G. Leal and A. del Ri'o in |CJLdR04) . Using the method 
of Remark 1351 our code can compute a fundamental domain for the group r(0). 
It computes first a maximal order C D O, and a fundamental domain for r(O') 
(having covolume 0.8889 . . . ). By coset enumeration, it finds that r(C') has index 9 
in the larger group, and computes a fundamental domain for the initial group r(C'). 
The overall computation takes 15 seconds. 

3.3. A larger example. Consider the ATR field F generated by an element t 
such that ti" + - \%t' - 211^ + 2U^ + 57^^ - 2t^ - ZZt^ - lOi + 1 = 0, having 
discriminant —546829505431 ~ —5.510^^. Let B be the quaternion algebra 

where a ^ \ [-2ht^ - 82t^ + 6lf + AOAt^ + 376t^ - 932^* - 718^^ + 590^^ + 368i - 
33) and b — ~1, which is ramified exactly at the real places of F. Let O be a 
maximal order in B. The group T{0) has covolume 1783.7 .... Our code computes 
a fundamental domain for this group in 23 hours and 39 minutes (85150 seconds). 
It spends 5.3% of the time for enumeration, 5.8% for the routine KeepSameCroup, 
87.7% for CheckPairing and 1.3% for computing the volume of the polyhedron. The 
fundamental domain has 5434 faces and 16252 edges. 

3.4. Efficiency of the algorithm. According to geometers, the parameter encod- 
ing the complexity of an arithmetic Kleinian group is the covolume. In practise it 
is simpler to vary the discriminant of the base field (and hence the degree) and the 
norm of the discriminant of the quaternion algebra. It seems hard to estimate the 
running time of the algorithm in terms of these parameters. First, we do not know 
any bound on the radii of the isometric spheres containing the faces of the Dirichlet 
domain, or of generators of the group, so we do not know how many elements we 
have to enumerate. Then, even if we have generators of the group, we do not know 
how long the normalized basis algorithm could run before terminating (see also 
Remark [?7|). 
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Figure 3.3. Running time of the algorithm 

We present some numerical data obtained in a family. Since the running time 
increases very quickly with the discriminant of the field, we fixed the base field and 
varied the discriminant of the algebra. The field we chose is the ATR cubic field of 
discriminant —23. We computed groups T{0) for every algebra with discriminant 
less than 10 000, and one algebra every ten with discriminant less than 15 000. 

Analysis of this data shows that the running time is approximately proportional 
to the square of the covolume, with a few exceptionnally slow computations. We 
explain this as follows: in almost all cases, the enumeration appears to take neg- 
ligible time, and the longest part is the computation of the fundamental domain 
itself; moreover the data (Figure I3.4|) seem to indicate that the number of faces is 
proportional to the covolume (we have such a lower bound since the volume of a 
hyperbolic tetrahedron is bounded by 3£(^)), and we know that our algorithm to 
compute the domain given the faces is quadratic. 
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Figure 3.4. Number of faces of the Dirichlet domains 
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